2. Create Mosaic of 2011 NDVI ImagesΒΆ
Written by Men Vuthy, 2021
Import modules
[1]:
import rasterio
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
import glob
import os
Create path of all ndvi image files
[2]:
# File and folder paths
file_path = "input/ndvi_2011"
# Make a search criteria to select the ndvi files
q = os.path.join(file_path, "20*.tif")
ndvi_fp = sorted(glob.glob(q)) # sorted files by name
print(ndvi_fp)
['input/ndvi_2011\\2010_01_01.tif', 'input/ndvi_2011\\2010_01_17.tif', 'input/ndvi_2011\\2010_02_02.tif', 'input/ndvi_2011\\2010_02_18.tif', 'input/ndvi_2011\\2010_03_06.tif', 'input/ndvi_2011\\2010_03_22.tif', 'input/ndvi_2011\\2010_04_07.tif', 'input/ndvi_2011\\2010_04_23.tif', 'input/ndvi_2011\\2010_05_09.tif', 'input/ndvi_2011\\2010_05_25.tif', 'input/ndvi_2011\\2010_06_10.tif', 'input/ndvi_2011\\2010_06_26.tif', 'input/ndvi_2011\\2010_07_12.tif', 'input/ndvi_2011\\2010_07_28.tif', 'input/ndvi_2011\\2010_08_13.tif', 'input/ndvi_2011\\2010_08_29.tif', 'input/ndvi_2011\\2010_09_14.tif', 'input/ndvi_2011\\2010_09_30.tif', 'input/ndvi_2011\\2010_10_16.tif', 'input/ndvi_2011\\2010_11_01.tif', 'input/ndvi_2011\\2010_11_17.tif', 'input/ndvi_2011\\2010_12_03.tif', 'input/ndvi_2011\\2010_12_19.tif', 'input/ndvi_2011\\2011_01_01.tif', 'input/ndvi_2011\\2011_01_17.tif', 'input/ndvi_2011\\2011_02_02.tif', 'input/ndvi_2011\\2011_02_18.tif', 'input/ndvi_2011\\2011_03_06.tif', 'input/ndvi_2011\\2011_03_22.tif', 'input/ndvi_2011\\2011_04_07.tif', 'input/ndvi_2011\\2011_04_23.tif', 'input/ndvi_2011\\2011_05_09.tif', 'input/ndvi_2011\\2011_05_25.tif', 'input/ndvi_2011\\2011_06_10.tif', 'input/ndvi_2011\\2011_06_26.tif', 'input/ndvi_2011\\2011_07_12.tif', 'input/ndvi_2011\\2011_07_28.tif', 'input/ndvi_2011\\2011_08_13.tif', 'input/ndvi_2011\\2011_08_29.tif', 'input/ndvi_2011\\2011_09_14.tif', 'input/ndvi_2011\\2011_09_30.tif', 'input/ndvi_2011\\2011_10_16.tif', 'input/ndvi_2011\\2011_11_01.tif', 'input/ndvi_2011\\2011_11_17.tif', 'input/ndvi_2011\\2011_12_03.tif', 'input/ndvi_2011\\2011_12_19.tif', 'input/ndvi_2011\\2012_01_01.tif', 'input/ndvi_2011\\2012_01_17.tif', 'input/ndvi_2011\\2012_02_02.tif', 'input/ndvi_2011\\2012_02_18.tif', 'input/ndvi_2011\\2012_03_05.tif', 'input/ndvi_2011\\2012_03_21.tif', 'input/ndvi_2011\\2012_04_06.tif', 'input/ndvi_2011\\2012_04_22.tif', 'input/ndvi_2011\\2012_05_08.tif', 'input/ndvi_2011\\2012_05_24.tif', 'input/ndvi_2011\\2012_06_09.tif', 'input/ndvi_2011\\2012_06_25.tif', 'input/ndvi_2011\\2012_07_11.tif', 'input/ndvi_2011\\2012_07_27.tif', 'input/ndvi_2011\\2012_08_12.tif', 'input/ndvi_2011\\2012_08_28.tif', 'input/ndvi_2011\\2012_09_13.tif', 'input/ndvi_2011\\2012_09_29.tif', 'input/ndvi_2011\\2012_10_15.tif', 'input/ndvi_2011\\2012_10_31.tif', 'input/ndvi_2011\\2012_11_16.tif', 'input/ndvi_2011\\2012_12_02.tif', 'input/ndvi_2011\\2012_12_18.tif']
List of ndvi image date and image array
[3]:
img_date = []
for file in ndvi_fp:
name = os.path.basename(file)
name_ver2 = name.replace("_","/")
name_ver3 = name_ver2.replace(".tif","")
img_date.append(name_ver3)
print(img_date)
['2010/01/01', '2010/01/17', '2010/02/02', '2010/02/18', '2010/03/06', '2010/03/22', '2010/04/07', '2010/04/23', '2010/05/09', '2010/05/25', '2010/06/10', '2010/06/26', '2010/07/12', '2010/07/28', '2010/08/13', '2010/08/29', '2010/09/14', '2010/09/30', '2010/10/16', '2010/11/01', '2010/11/17', '2010/12/03', '2010/12/19', '2011/01/01', '2011/01/17', '2011/02/02', '2011/02/18', '2011/03/06', '2011/03/22', '2011/04/07', '2011/04/23', '2011/05/09', '2011/05/25', '2011/06/10', '2011/06/26', '2011/07/12', '2011/07/28', '2011/08/13', '2011/08/29', '2011/09/14', '2011/09/30', '2011/10/16', '2011/11/01', '2011/11/17', '2011/12/03', '2011/12/19', '2012/01/01', '2012/01/17', '2012/02/02', '2012/02/18', '2012/03/05', '2012/03/21', '2012/04/06', '2012/04/22', '2012/05/08', '2012/05/24', '2012/06/09', '2012/06/25', '2012/07/11', '2012/07/27', '2012/08/12', '2012/08/28', '2012/09/13', '2012/09/29', '2012/10/15', '2012/10/31', '2012/11/16', '2012/12/02', '2012/12/18']
[4]:
# List for the source files
img_array = []
# Iterate over raster files and read them
for path in ndvi_fp:
image = rasterio.open(path)
image = image.read(1)
img_array.append(image)
Normalize ndvi images
[5]:
def cal_index(img):
# By default numpy will complain about dividing with zero values.
# We need to change that behaviour because we have a lot of 0 values in our data.
np.seterr(divide='ignore', invalid='ignore')
# Convert ndvi scale index to -1, 1
ndvi_idx = img/10000
return ndvi_idx
[6]:
# Normalize all images
normalized_image = []
for array in img_array:
normalized_image.append(cal_index(array))
NDVI_Images = np.array(normalized_image)
# Check Shape
NDVI_Images.shape
[6]:
(69, 521, 753)
Save as Raster
[7]:
# Add one image for projection and shape reference
raster = rasterio.open("input/ndvi_2011/2011_07_28.tif")
[8]:
# Data dir
data_dir = "output/1/mosaic_img/"
# Output raster
out_tif = os.path.join(data_dir, "NDVI_2011.tif")
# Copy the metadata
out_meta = raster.meta.copy()
out_meta
# Update the metadata
out_meta.update({'driver': 'GTiff',
'dtype': 'float32',
'nodata': None,
'width': raster.shape[1],
'height': raster.shape[0],
'count': NDVI_Images.shape[0],
'crs': raster.crs,
'transform': raster.transform
})
with rasterio.open(out_tif, "w", **out_meta) as dest:
dest.write(NDVI_Images.astype(np.float32))
Save ndvi plot into PNG files
[9]:
for i in range(len(img_date)):
plt.rcParams["figure.figsize"] = (10,10)
plt.imshow(normalized_image[i], cmap='PiYG', vmin = -1, vmax = 1)
plt.colorbar(fraction=0.03, pad=0.04)
outfp = 'output/1/png_img/'
plt.title(img_date[i])
plt.savefig(outfp+'{0:.0f}.png'.format(i), dpi = 300)
plt.close();
Create GIF files to see timelapse of ndvi from 2010-2012
[10]:
import natsort
from PIL import Image
fp_in = "output/1/png_img/*.png"
fp_out = "output/1/gif_img/ndvi_change_2010_2012.gif"
png_unsorted = glob.glob(os.path.join(fp_in))
png_img = natsort.natsorted(png_unsorted,reverse=False) # sorted files by name
img, *imgs = [Image.open(f) for f in png_img]
img.save(fp=fp_out, format='GIF', append_images=imgs,
save_all=True, duration=200, loop=0)
[11]:
from IPython.display import Image
Image(filename="output/1/gif_img/ndvi_change_2010_2012.gif.png")
[11]:
